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Methods of dynamical systems have been used to study homogeneous and isotropic cosmological 
models with a varying speed of light (VSL). We propose two methods of reduction of dynamics to 
the form of planar Hamiltonian dynamical systems for models with a time dependent equation of 
state. The solutions are analyzed on two-dimensional phase space in the variables (a;, x) where x 
is a function of a scale factor a. Then we show how the horizon problem may be solved on some 
evolutional paths. It is shown that the models with negative curvature overcome the horizon and 
flatness problems. The presented method of reduction can be adopted to the analysis of dynamics of 
the universe with the general form of the equation of state p = y(a)e. This is demonstrated using as 
' an example the dynamics of VSL models filled with a non-interacting fluid. We demonstrate a new 

type of evolution near the initial singularity caused by a varying speed of light. The singularity-free 
oscillating universes are also admitted for positive cosmological constant. We consider a quantum 
VSL FRW closed model with radiation and show that the highest tunnelling rate occurs for a 
constant velocity of light if c(a) oc a n and — 1 < n < 0. It is also proved that the considered class of 
models is structurally unstable for the case of n < 0. 
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PACS numbers: 98.80.Cq, 95.30.Sf 



I. INTRODUCTION 



Although the standard cosmological model is usually believed to be a correct picture of our universe [1| it has 
still some difficulties among which the flatness and horizon problems are most widely known. The existence of the 
large-scale structure in the universe extending to the limit of the deepest surveys is another mystery. Its very presence 
O implies the appearance of some seeds for this structure in the early universe. In the standard big-bang scenario they 
!— i . should be built in, which is rather an undesired feature of the theory. Therefore, so great attention has been paid 
to inflationary universe models which, albeit invoking an exotic (if not hypothetical) physics, were able to provide at 
, least some hope for a consistent explanation of both flatness and horizon problems as well as the origin of seeds for the 
T ■ large-scale structure. The results of the early universe physics lead us to expect the occurrence of phase transitions 
r> \ when the universe was young, hot, and dense. 

The varying speed of light (VSL) cosmology, seen as an alternative to the inflation theory, was proposed by Moffat 
0, who conjectured that a spontaneous breaking of the local Lorentz invariance and diffeomorphism invariance 
associated with a first order phase-transition can lead to the variation of the speed of light in the early universe. This 
idea was revived by Albrecht and Magueijo Q and was given further consideration by Barrow 0-0 Barrow showed 
that the conception of VSL can lead to the solution of flatness, horizon, and monopole problems if the speed of light 
falls at an appropriate rate. The dynamics of VSL was widely studied in theoretical as well as empirical contexts 
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The main motivation for the study of VSL models is to seek explanations for some unusual properties of the universe 
and to overcome some of the shortcomings of inflation scenario ^3 ■ I n particular, there is an empirical evidence of 
varying with time the fine structure constant in the context of the consistency of quasar absorption spectra |18| . 
Moreover, unlike the inflation the VSL theory provides a solution of the cosmological constant problem. However, it 
cannot solve the isotropy problem. It is also interesting to evaluate the power of this model to explain the acceleration 
problem [313. 

Of course, the VSL model as well as other models discussed in the literature have an ad hoc element (variable c) 
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not yet firmly founded within any existing physical theories. This feature does not seems to be so exotic to discard 
these models from discussion in the scientific community. Some brilliant arguments justifying this approach are given 
by Albrecht and Magueijo Q. 

The varying speed of light models which provide decent fits to the real universe are characterized by the speed of 
light or gravitational coupling which varies with time in the very early universe but is nearly constant today. Because 
there are stringent bounds on how fast these constants can vary with time after the first few seconds, the models 
which dynamics we study in the paper are only relevant in the very early universe. This should be made very clear 
by using the phase space approach and its tools for classi fying the qualitative types of solutions [2lt . 

The present paper is a continuation of previous papers [23, UM 011 * ne dynamics of VSL cosmology. We introduce 
a simple framework which allows to study the dynamics of VSL models in a general way, independent of any specific 
assumption about the equation of state, or the behaviour of the scale factor a(t) near the spatial singularity. We 
formulate a VSL Friedmann-Robertson- Walker cosmology as a two-dimensional dynamical system and we discuss its 
properties on phase portraits where trajectories represent all solutions for all physically admissible initial conditions. 
The methods of dynamical systems allow to indicate how the existence of certain desired physical effects depends on 
the choice of initial conditions and to analyse how these initial conditions determining the corresponding solutions 
are distributed in the phase space. 

Our main goal was to perform a global analysis of the dynamics of VSL cosmological models. We avoid the as- 
sumption of power type evolution in the VSL models, which are represented by critical points (singular solutions) in 
the phase space. We analyse the dynamics of the models on the phase plane and discuss how different trajectories 
representing non-singular solutions can solve cosmological puzzles. We conclude that models with the negative curva- 
ture and positive cosmological constant are preferred (in the sense that they have the largest set of initial conditions 
leading to the solution of flatness and horizon problems) . 

On the other hand we present two arguments which distinguish the FRW models with constant velocity of light. 
The theoretical one is that the VSL FRW models are structurally unstable if c(a) cx a n and n < contrary to classical 
FRW models. The quantum mechanical one is that if a closed universe was born from a quantum fluctuation via the 
quantum tunnelling process then the most probable universe is that with c = const. In this interval the potential 
function preserves its classical character and the universe tunnels from a zero size. 

The dynamics of considered cosmological models is reduced to the dynamics of a unit mass particle in one- 
dimensional potential. Then different physical properties like flatness, horizon and cosmological contant problems 
can be formulated in terms of the diagram of potential function of the system. 



II. THE METHOD OF THE DYNAMICAL SYSTEM STABILITY 



First of all, equations describing a cosmological model should be reduced to the form of a dynamical system 

x-i — fi (x\ , . . . , x n ) , z 1, . . . , n 
at 

in such a way that the solution with a static microspace (or other solutions of interest) is a critical point of the system 
(a;*,.. . ,<), i.e. for every i, /i(xj,. . . , a;*) = (i = 1, . . . , n). 

If a critical point is non-degenerate, i.e., at this point all real parts of the eigenvalues ReA^ of the linearization 
matrix 

A*. = M 

3 dx-j 

J x=x* 

do not vanish, then there is a one-to-one continuous mapping of a neighbourhood of this point which transforms 
trajectories of the original system into trajectories of the linearized system 




In this sense, the qualitative behaviour of the original system is equivalent to the behaviour of its linearized part. If 
(£<) • ■ ■ i£f) are eigenvectors of the linearization matrix A*, the solution of the linearized system has, in general, the 
following form 

n 

Xi (t)-x* =Re^C fc ^e Afc * 
fc=i 
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where Ck are constants. A non-degenerate critical point is called an attracting point if, for all eigenvalues, Re A; < 0. 
In this case, all trajectories from a neighbourhood of this point go to it when t — ► oo. A non-degenerate critical point 
is called an repulsing point if, for all eigenvalues, Re A; > 0. In this case, all trajectories from a neighbourhood of 
this point go to it when t — > — oo. A non-degenerate critical point is said to be an unstable saddle point and if a 
dynamical system has, at (x* ,...,£*), d eigenvalues with positive real parts and n — d eigenvalues with negative real 
parts (d = 1, . . . , n — 1). 

When investigating the stability of solutions with a static microspace the following theorem proves to be of special 
interest. If x* is a non-degenerate critical point and if the dynamical system has, at x* , d eigenvalues Ai, . . . , Ad with 
negative real parts, then there exists (locally) an invariant d-dimensional stable manifold Wf t , on which all trajectories 
of the system go to x* as t — > oo. (A manifold M is said to be an invariant manifold of a system if every trajectory 
passing through a non-degenerate point of M lies entirely in M (for — oo < t < +oo). For every such solution, there 
exists the asymptotic 

f " 1 

lim r 1 In < ^2[xj(t) - x*} 2 \ = a t (1) 



for a certain i. Similarly, if at the critical point x* the system has k eigenvalues with positive real parts then there 
exists an invariant fc-dimensional unstable manifold W„ nst on which all trajectories go away from the critical point 

m 

From the above theorem follows that, for a saddle point, there are two invariant manifolds Wf t and W™ ns f containing 
this point and filled with trajectories (separatrices) going to and going away from the critical point. All other 
trajectories (not contained in Wf t or in W^~ s ^) do not meet the critical point in question. 

For the complete construction of phase portraits in a plane it is necessary to know how trajectories of dynamical 
system behave at infinity. Let us take as an example the two-dimensional system 

x = P{x,y) (2) 
y = Q(x,y). (3) 

In the case of polynomial right-hand sides one usually introduces projective coordinates, eg. z = 1/x, u = y/x 
or v = l/y, w = x/y. Two maps (z,u) and (v,w) are equivalent if and only if u ^ and v ^ 0. Infinitely 
distant points of [x, y)-plane correspond to a circle S 1 which can be covered by two lines z = 0, — oo < u < oo and 
w = 0, — oo < v < oo. The original system in the projective coordinates (z, u) and after the time reparametrization 
r — * T\ : dr\ = xdr assumes the form 

z = zP*{z,u) (4) 
u = Q*(z,u) — uP*(z,u) (5) 

where 

P*(z,u) = z 2 P(l/z,u/z) 
Q*{z,u) = z 2 Q(l/z 7 u/z) 

and dot denotes differentiation with respect to new time t\ . 

On a similar way, in the projective coordinates (v, w) and in the new time T2 : dTi — ydr we obtain 

v = ~vQ* (v, w) (6) 
w = P* (v,w) — wQ* (v,w) (7) 



where 

P*(v,w) = v 2 P(l/v,w/v) 
Q*{v,w)=v 2 Q(l/v,w/v) 

and dot denotes differentiation with respect to time t%- 

The idea of structural stability was introduced by Andronov and Pontryagin [25| . A dynamical system S is said to 
be structurally stable if there exist dynamical systems in the space of all dynamical systems which are close, in the 
metric sense, to S or are topologically equivalent to S. Instead of finding and analyzing an individual solution of a 
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model, a space of all possible solutions is investigated. A given physical property is believed to be 'realistic' if it can 
be attributed to large subsets of models within a space of all possible solutions or if it possesses a certain stability i.e., 
if it is also shared by a slightly perturbed model. There is a well established opinion among specialists that realistic 
models should be structurally stable. What does the structural stability mean in physics? The problem is in principle 
open in higher than 2-dimensional case where according to Smale there are large subsets of structurally unstable 
systems in the space of all dynamical systems [2f|. For 2-dimensional dynamical systems, as in the considered case, 
the Peixoto's theorem says that structurally stable dynamical systems on compact manifolds form open and dense 
subsets in the space of all dynamical systems on the plane. Therefore, it is reasonable to require the model of a real 
2-dimensional problem to be structurally stable. 

In our further considerations we will investigate the dynamics of dynamical systems in a finite domain of phase 
space as well as at infinity. At this point we would like to recommend t he p resentation of the actual state of the art 
in the field of application of the dynamical systems to general relativity |22| • 



III. BASIC EQUATIONS OF THE THEORY 

Albrecht and Magueijo Q and Barrow |f| set up a useful framework to discuss the VSL models assuming that the 
time variable c should not introduce changes in the curvature terms of the gravitational field equations and that the 
Einstein equations must hold. Because varying c breaks the Lorentz invariance the VSL cosmology requires a specific 
reference frame (including a specific choice of a time coordinate) in which changes in the field equations are minimal 
and one postulates it to coincide with the cosmological comoving frame. 

In the case of the VSL version of the FRW models (with A = 0) the scale factor obeys the following dynamical 
equations 

d\ 2 8nG{t)p Kc 2 {t) 

( 8 ) 



a) 3 a 2 {t) 



a 



Equation Q is called the Raychaudhuri equation and from the above system one can build a generalized conservation 
equation 

. h ( p \ G ZKc 2 c 



c 2 (t)J ' G 8irGa 2 

in which time dependence of fundamental constants was taken into account explicitly. Alternatively one can think of 
the Raychaudhuri equation together with the generalized conservation equation as of a fundamental system to which 
equation J5J) is a first integral. 

The fundamental difficulty concerning system JSJ-^UJ) is that it is a non-autonomous system with unknown functions 
G(t) and c(t). In order to be specific in the further analysis we adopt Barrow's power-law Ansatz 

G(t) = G a(t) q , c(t) = c a(t) n . (11) 

Moreover, we assume the hydrodynamical energy- momentum tensor with the equation of state for the non- interacting 
multifluid 



e = 1 (a)pc 2 (12) 



where pi — pioa~ 3 ( li+1 \ and energy density e = pc 2 . 

In the special case of the matter and radiation mixture, the factor 7(a) depending on the scale factor a takes the 
form 

^( a ) = I — TT> P = °+l e ^ e = e m + e r (13) 
3 aa + 1 3 

where a = p m o/prO- 

If we substitute I = 1 into i|12|) then we obtain models filled with single matter and with the equation of state 
p = jpc 2 . Generally, the equations of state for non-interacting fluids with pressure p = '^2 i ^fiPiC 2 , the equation of 
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state assumes the form p — r y(a)pc 2 where the factor 7 can be parametrized by the scale factor. This fact is crucial 
for the reduction procedure. 

Power-law Ansatz turns the field equations back into an autonomous system. Now we can think about 
extensions of the baseline equations. One can in a straightforward way include the cosmological constant A by 
introducing pressure p\ and energy density p\ 

PA = -pAC 2 {t) (14) 
Ac 2 (t) , , 

" = **m- (15) 

System with the cosmological constant can be cast into form 

8nG(t)p Kc 2 (t) Ac 2 (t) 



a J 3 a 2 (t) 3 

a _ AnG{t) ( 3p \ Ac 2 (t) 



a 3 V c 2 (t) J 3 



(16) 
(17) 



'at v 9, ,\ G 3Kc 2 a , . 

p + 3- [p+^c 2 (t )=_p_ + _ --. 18 

Equation (|18|l is easy to solve only for the case oic — const. Therefore, in our consideration of the general formulation 
of dynamics we cannot use an explicit form of a solution of equation (|18|) . To avoid this difficulty we consider a special 
procedure of reduction. 



IV. REDUCTION TO A PLANAR HAMILTONIAN SYSTEM 



A. The general solution of a dynamical problem 



To construct a dynamical system we assume the form of the equation of state p = j(a)pc 2 and calculate the density 
p using both equation l(T3|l and (JTSJ). For simplicity we focus our attention on the case of E = corresponding to the 
VSL model with matter in the multifluid form. Then we obtain from equation (|16fl 



and from equation Ijl7(l 



8irGp 



8nGp 



Kc 2 (t) 



Ac2_ 

1T 



1 



2a 2Ac 2 



3 l + 37(a) V a 3 

By adding the sides of above equations we obtain a second order nonlinear equation with respect to variable a 



where 



n(a) 



a + %jj{a)a 2 + n{a) = 



1 + 37(a) 
2a 

^(1 + 37(a)) -|a(l+ 7(a)) 



c 2 (a) 



and only the term n(a) depends on the cosmological constant. 

Equation (|21|l can be rewritten as an autonomous dynamical system 



a = p 

p = —i/j(a)p 2 — n(a). 



(19) 



(20) 



(21) 



(22) 
(23) 
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To apply the dynamical system theory, it is useful to reduce system ()22|l - (|23|l to the form with polynomial right-hand 
sides 



, da 
a = — — pa 

07) 

P = ^--i(3 7 (a) + l)p 2 -0(a) 



(24) 
(25) 



where 4>{a) — an{a) and t — * rj: — — dr). The solution of equations l|24|l -H25 |) represents a phase curve in the phase 
space (a,p) = (a, a). 

The solution of equation (|21|l may be given after substitution a = p(a) or equivalently by taking the quotient of 
equations (12211 and (|23|> . Then we obtain 



p^- + 4>{a)p 2 + k(o) = 0. 
da 



(26) 



Equation (1201) takes the form of the Bernoulli equation and after standard substitution u(a) — p 2 , u' = du/da we 
obtain the non-autonomous system 

u + 2ij)(a)u + 2k(o) = 0. (27) 
Finally, the solution of equation l|21f> . passing through the point (ao,po(oo))> can be given in the following form 



u(a) = cxp —2 / ip(a)da 

J ao 



Poiflo) ^ I 2k(o) exp ( 2 j ip(a)da ) da 

ao 



p 2 {a) = exp I - / 

\ J a, 



37(a) + 1 



da \ Po(ao) - 2 



K(l + 37(a)) 
2a 



Aa( 7 (a) + 1) 



exp 



a 3 7 (a) + 1 



da da 



To consider the case of a mixture of matter and radiation, we substitute the special form of 7(a) from formula 1|13|) 
and then we obtain 



ag(aa + 1) 



a z (aa + 1)1 



aa + 1 



-da 



K(l + 37(a)) 



2a 



Aa 1 



1 



3(aa + 1) 



c 2 (a) 



and the general solution of equation (|27() has the form 



u(a) = exp I —2 / ip(a)da 



C — 2 / k(o) exp 2 / ip(a)da I da 



(28) 



(29) 



It means that the following expression can be treated as a first integral of system (|24|l - l|25() . It is characteristic for 
dynamical systems of general relativity and cosmology that a first integral can be used in constructing a Hamiltonian 
function. A first integral can be represented as algebraic curves in the phase space. These algebraic curves are given 
by 

(30) 



p (a) exp ( 2 / ip(a)da I + 2 / n(a) exp I 2 / if) (a) da ) da = C. 

J ao J J ao \ J ao 



And then in the considered case we obtain 



■■^— + 2 / K(a)^— 
aa + 1 /.„ aa + 1 



da = C 



(31) 



Now if we introduce a new variable, x, such that 



1 ada 
—=dx = ■ 
V2 y/aa + 1 



the above relation can be written as 



+ V(x) — V(a ) = const 
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where 

plays the role of the potential with a(x) : -^=x = -^\Jaa + l(aa — 2). 

This procedure successfully works for any function ip(a). It is sufficient to replace the expression a ° +1 in the 
considered case by exp(2 f \jj{a)da) = 2 . 

In the special case of a = we obtain that p = and we can use the standard formalism considered in |2Sj . 

Now we can see from formula Ij28f) that an algebraic curves on which lie trajectories of the system take the com- 
plicated form. Therefore it is useful to visualize them in the phase space. Using the form of the above first integral 
we can also classify all possible solutions by considering a limiting curve x — and derive the relation A(a) as was 
presented in the classical case pflj . 

There are two kinds of the different methods of reducing equation Ij21|l to the form of Newtonian equation of motion 
in a one-dimensional configuration space. First, after introducing the new rescaled variable a — > x, we obtain dynamics 
in the form x = — Second, after introducing the new time variable, say r(t) defined in such a way that the term 

ip(a)a 2 can be dropped in this parametrization and we obtain dynamics in the form x" = 4j5 = — 

Let us note that the function V(a(x)) plays the role of the potential for a particle which the position is given by x 
and motion described by 

.. _ dV(a(x)) 
dx 

Szydlowski et al. 28] showed that the choice of new variables x = x(a) allows to reduce the dynamics of classical FRW 
models with matter or radiation to a one-dimensional Newtonian equation of motion. It is also interesting that there 
is a similar possibility of reducing system l|21|l to the Hamiltonian form for any equation of state 7 = 7(a), for example 
for any mixture of non-interacting fluids. To perform this let us consider the general nonlinear reparametrization of 
variable a such that 

x = a D ^ a » = a D ^ (33) 

It can be shown that equation l|21|l can be reduced to the form of a Newtonian equation of motion, 'x = —dV(x)/dx, 
if the multiplicative coefficient appearing in d 2 vanishes. This condition gives us the following 

D aa (hxa) + — ( J D a lna + — = V>(a) ( D a In a + -) (34) 
a a z \ a J \ a J 

where D a = dD/da. In the special case of 7(a) = const (D a = 0), classical results can be recovered |28|: we have 
D = 2 for pure radiation and D = 3/2 for dust and D = |(1 + 7) for perfect fluid with p — 7/9. If 7 = 7(a) is any 
function of a then D(a) is a solution of the above equation. After simple substitution 

D a (\na) + ^- = z(a) (35) 
a 

we obtain that z(a) is a solution of the equation 

*+ — = #*)• (36) 
z da 



In this new variable 

x = —K,(a{x))z{a{x))x 



dV 

dx 



and the dynamics is reduced to the case of a nonlinear oscillator with 'spring-like tension' k(x) — K,(a)z(a). 

The information about the equation of state is hidden in the function ip(a) and after finding the solution z{a) for a 
specific form of the equation of state it should be easy to find D(a) from equation (|35|l . It can be easily shown that 
the corresponding equation determining z(a) is the Bernoulli equation for which the solution is 

z(a)= 7^k = ^ ln (/ Ha)da 
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where </>(a) = exp(J %jj(a)da). If we put z(a) into (|35[1 we can find that 

, . In f a d>(a)da f a , % 

D(a)= In a =1 °g«7 ^°) da - 

For the case of mixture of radiation and dust it has the simple form 



(aa+l)(aa-2) 

z ( a ) = { #■ for a = oo 



for a 7^ 
for a = oc 
for a = 

where </>(a) = ^ a ° +1 and f a (j)(a)da = -^\\/aa + l(aa — 2)|. 

Therefore we obtain for the case of non-interacting matter and radiation 





— 3^ L J for a 

^(a) = { 2 for a = 

for a = oo. 

It can be proved that in a general situation we have the following relation 

<j)(a)da = x(a) = / Vae^ J ° ^ da ') da (37) 
and 

f a lo g^ 2 ^^^ — ) fora^O 
»(o) = < a 2 for a = 

[a 3 / 2 for a = oo 

i.e., for any fluid (or its mixture) which satisfies the equation of state for so-called 'quintessence' matter p — r y(a)pc 2 (a) 
we can always find the corresponding D(a). 

Due to 13 7|) the equation of motion can be rewritten to the simplest form 

x = —K(a(x))<f>(a(x)), V(x) = / K(a(x))4>{a{x))dx. (38) 

Let us note that there is also a possibility to generalize such a result to the case of non- vanishing shear in B(I) or 
B(V) models when a oc x z / 2 . 

In the second approach it is useful to reparametrize the time variable t in such a way that 

t -> t : dt = cj)(a(T))dT (39) 

where ip(a) is a yet-to-be-determined function which should be chosen in such a way that term ip(a)a 2 is absent in 
1)21(1. We can do that provided that <f> satisfies condition 



(j) ee exp / tp(a)da. (40) 

Then equation (|21|l assumes the form similar to the equation for the motion of a non-relativistic particle in the external 
field with the potential V(a), namely 

^°=-^ (41) 

and 

V(a) = f K(a)cf> 2 {a)da (42) 
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where 



0(a) 
k(o) 



1 + 37(a) 

2a 

£(1 + 37(a)) -^a(l+ 7(a)) 



c 2 (a). 



Therefore for such a system the Hamiltonian takes the form 



H{a,p) 



+ V(a) = E = const 



(43) 



where the correspondence to the vacuum case is reached after putting E — and 7(a) = 0. 

The advantage of this procedure of reduction is its simplicity. The new time variable r is a monotonic function of 
Newtonian time t and motion is represented in the form of one-dimensional Hamiltonian system with the potential 



V(a) = / K(a)4> 2 (a)da = 4 



^(1 + 37(a)) -|a(l+ 7(a)) 



a 2n <t> 2 {a)da. 



(44) 



By comparison potential \Y'>'1\\ and potential (|44|l we can observe that both procedures give rise to the same form 
of the potential function as a function of a. Whereas the second approach seems to be simpler, the first one has an 
advantage that it allows to discuss dynamics in the origin time. In both cases we do not explicit integrate continuity 
equation l|18|) which gives us the relation p(a). The effect of matter content is included in p{a) and energy constant 
E. 

For the special case of matter content in the form of non-interacting matter and radiation we obtain 



K4 [ a (aa + 2)a 2 "+ 1 
W 2 J (aa + iy 



da 



Ac 2 f a (3aa + 4)a 2 "+ 3 
(aa+ l) 2 



da. 



(45) 



where integrals in the above form of the potential can be given explicitly. 
In the special case of A = n = we have 



V(a) 



Ka 



2(aa + l) 

Therefore the dynamics is given by the Hamiltonian equations 

, &H 

dp 



P = 



on 

da 



(46) 

(47) 
(48) 



which constitute the two-dimensional dynamical system. Entire evolution is represented by an evolutional path on 
the plane (a,p). The domain of acceleration a(t) > is determined by the condition 



The above condition can be rewritten as 



a" - (a') 2 V'(a) > 0. 



r)V 

~ - (a')V(a) > 0. 
oa 



(49) 



(50) 



Let us note that equation (|50|l can be formulated equivalently as a following condition in some domain of the config- 
uration space {a: a > 0} 



~-2(B-V(aM«)>0 
oa 



or 



where 



1 dV 

■0(a) da 



2V{a) < -2E 



(51) 



(52) 



dV 
da 



-K{a)4> 2 {a) and V(a) = 1+ ^ (a) . 
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V. EVOLUTION OF THE VSL DYNAMICAL SYSTEM ON PHASE DIAGRAMS 



A. Background 

In the further qualitative analysis of dynamical system H22(l - (|23[) we consider the matter as the mixture of radiation 
and dust. Then system l|22|l - (|23[) takes the form of the autonomous system with rational right-hand sides 

a = p (53) 

aa + 2 2 / K(aa + 2) 3aa + 4 \ 2 , 

-p — I —— -Aa)c 2 (a) (54) 



2a(aa + l) \2a(aa + l) 6(aa + 1) 

where c 2 (a) = c§a 2 ™, n < 0, a = e m o/e r o (for regularization of the system in a = it is useful as in (|24|I - H25|) to 
introduce time r\ : ^ = drf) . Of course the above system possesses the first integral in form (|28|) . 

In the finite domain, system H53|) - (|54[) has at most one critical point which corresponds to an extremum of the 
function V(a) : ^| a _ ar = 0, po(ao) = 0. The stability of this point is determined from the convectivity of a diagram 
of potential function V(a). There are two limit cases corresponding to the equation of state: of dust (a — » oo) and of 
pure radiation (a = 0). From system (|53|) - (|54|l in the latter case we obtain the VSL system with pure radiation 

a = p (55) 
V -(-- l^a) c 2 {a). (56) 



a \ a 3 

In the above system instead of a in the spirit of the first approach we introduce a new variable x = a 2 and we 
obtain the system in a simpler form 

x = y (57) 
y=-2^K-^Ax^jc 2 (a(x)) (58) 

where c 2 (a(x)) = CqO 2 " = Cqx". The phase portraits of system (|57|l - (|58|l for n = —2 are presented on Fig. ^ 
In the projective coordinates (z,u) the above system takes the form 

z = —zu (59) 
u = -2 \Kz - ^A^j c 2 z- n - u 2 (60) 

This form of the system is useful in analysis of behaviour of trajectories at infinity, because z = corresponds a circle 
at infinity x = oo which bounds the phase plane. The phase portraits of system (|59() - (|60|l for n = —2 are presented 
on Fig. 

In this case first integral (|28|l takes the form 

2 - N_2, 



— + 2 J yK- -Axj c 2 (x)dx = C = const. (61) 

Let us note that in the special cases of n = —1 and n = —2 the potential takes the particular form V{x) = 
2{K\nx - fAz), V(x) = 2(—| - fAlnx). 

It is clear that first integral (j61(l is in fact the integral of energy because system (|57|I - H58|I is a Hamiltonian dynamical 
system with the Hamiltonian 

■n 2 ( Kr n+1 2 T ™+ 2 \ 

*<»*> = IT + 2 (— - 5 A ^) s c - »' I8t > 

Now the integral of energy can be used in the classification of all possible evolution modulo their quantitative (i.e., in 
accuracy to differential type) properties. 

In the qualitative classification, for any case of 7(a), first integral l|28|l may be useful as in the method of the 
classification previously used in [2|| ■ It assumes the form 

V 2 

y — + V(a(x)) = V(a ) (62) 
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K = -l, A=l 



K = -l, A = 








FIG. 1: The phase portrait for system JHEJ for n = -2. For the cases (K = -1, A = 1), (A" = -1, A = 0), (if = 0, A = 1) 
the qualitative structure of the phase space is the same. In the finite domains of phase space there are no critical points. Typical 
trajectories represent the solution starting from singularity-free stage x = oo at t = 0, then reach the stage x = xo,y — 
and goes to infinity. In this three cases the trajectories pass through the acceleration region (the shaded region). The largest 
acceleration region is for the (K — — 1,A = 1). These models accelerate for finite interval of time, and this acceleration 
happens to the models without the cosmological constant. The closed models for (K = 1, A = 0) are the typical oscillating 
models. There are no acceleration for (K = 1, A = 0), (K = 0, A = 0). In the case of (K — 1,A = 1) the saddle appears 
and we have two additional types of evolution. The Lemaitre model with quasi-static stage of evolution and models realized in 
the aforementioned opened and flat models. The acceleration region is in the middle of quasi-static phase of evolution of the 
Lemaitre. 



where 



V(a) 



a(x) 



1 



X(l + 3 7 (q)) 
2a 



Aa 1 



3(aa + 1) 



c 2 (a)da, 



V(a ) = const > and x = ^0-t 2 (t — 1), t = y/aa + 1 and x = a 2 for a = 0. 

From equation l|62(l . after imposing the condition p = 0, we can calculate A from the expression on function A(a) 
which constitutes a boundary of a domain of configuration space admissible for motion 



p 2 > V(a) - V(a ) < 
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FIG. 2: The phase portrait for system l|59(l - l)6()ft for n = —2 in the projective coordinates z = 1/a;, u = y/x. All critical 
points at infinity corresponds to zq — 0. There are two types of critical points (zo,uo) = (0,0) or (zo,uo) = (2A/3if, 0). 
For all points tr A = (where A is a linearization matrix of the considered system), i.e., at least one of the eigenvalues is 
zero. If (A > 0, K = 1) we have a saddle point at (20, 0). The presence of degenerate points at infinity (0,0) that the 
system is structurally unstable in contrast to the models with positive cosmological constant and constant velocity of light 
where (0, w 4A/3) represents a stable node (the de Sitter stage). For n > there is no critical points at (0, 0) and models are 
structurally stable. 



and 



A(o) > 



I a ^$KcHa)da-V(a ) 



(63) 



By consideration the boundary dD: {(a, A): A(a) = 0} and construction of the levels A = const we obtain the 
qualitative classification of all possible trajectories in the space (A, a) or (A, x). 
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B. How to interpret the acceleration of scale factors and the absence of particle horizon 

A great advantage of the phase-space dynamical description is the ability to discuss the distribution of models with 
given properties. In the other words one can imagine an ansamble of models starting from different initial conditions 
and ask how a given property is distributed in the ansamble. Now we formulate sufficient conditions for solving the 
flatness and horizon problems in terms of phase-space relations. Let us recall that the flatness problem is solved 
whenever the scale factor's acceleration is positive 

a(t) > 0. 

This condition is fulfilled in the subspace 2?fl a t of the phase space 

aa + 2 2 / K(aa + 2) iaa + 4 \ 2 



2a(aa+l) \2a(aa + 1) o(aa + 1) / J 

It means that trajectories representing the histories of VSL universes undergo an accelerated expansion while staying 
in region 2?fl a t- One can restate relation l|64|) using the Hamiltonian constraint p 2 = 2(V( a o) — V( a ))- It is easy to 
see that the respective condition expressed purely in terms of configuration space, reads 

-r, f era + 2 .... , ... .. fK(aa + 2) iaa + A . \ 2 . . 1 . c . 

^accei = \a: - —AV{a ) - V(a)) - — — 4- - — Aa c 2 a >0 . 65 

^ a{aa + 1) \2a(aa+l) o(aa+ 1) / J 

Analogous criteria of acceleration if dynamics is covered by equations (|47|l - l|48|l are given by (|50|) in phase space and 
(|52H in the configuration space. 

Another interesting question concerns the horizon problems. It is easy to prove the following criterion of avoiding 
the horizon problem. 

Theorem 1 The FRW cosmological model does not have an event horizon near the singularity if a(t)c~ l (t) tends to 
a constant while a(t) tends to zero. 

Proof. When all events whose coordinates at past time are located beyond some distance dn they can never commu- 
nicate with the observer at coordinate r = in the Robertson- Walker metric we can define the distance dn as past 
event horizon distance. It is given by 

, . . , , f to dt'M) 
d H (t) = a(t) —±±=a(t)I. 
Jo a{?) 

Whenever I diverges as t — ► there is no past event horizon in the space time geometry. On the other hand, when / 
converges the space time exhibits a past horizon 



to 



dtc{t) _ f to a n da _ f to 1 da 



Let c 1 (t)d < A then 



da 1 n N 
/ — = — (lna + oo) 
Jo a A 



On the other hand when a < A is bounded then 

<- to dtc{i) _ f to a n da _ f to n _ x da 



a(t) 



and I > i/ ao a" 1 ' 1 da or 



s ; t0 dtcjt) > a* 
a(t) ~ n 



Therefore / diverges as a — > if n < and a < A □. 
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The above criterion can be reformulated in the language of the phase space in the form 



a — > and c (t)a(t) — > const, i.e. a 



-2)1 



(const)^ 



For example for the radiation case a = y/x and then the past horizon is eliminated if only x 2 x~ 2n ~ 1 — > const as x — > 0. 
After substituting first integral (|61|l we obtain that there is such n that the horizon disappears if only if n < — 1. Let 
us note that the above proof is based on the Hamiltonian constraint and is independent of any specific assumption 
about an equation of state or a(t ) near the singuarity. If we assume the power-law behaviour of a(t) then Barrow's 
result can be simply achieved [l7|, namely a(t) oc £ 2 / 3 (t _1 ) anc j if 2 > 2n + 3(7 + 1). 

C. The potential function for mixture of dust and radiation 

The boundary of the domain in the configuration space (x-space or a-space) admissible for motion is determined 
by expression for 



V(a(x)) = 1(4 



J (aa+l) 2 °J 



a{x) (3aa + 4)a 3+2n 



,-,,-> -■ <>, .,, , , a - d a = V{a ) > 0. 

{aa+\y J 3(cra+l) z 



The equation V{a(x)) — V{ao(xo)) — can be represented in the space (A, a) or more correctly in the space (A, x) as 
a boundary curve for the classification. Instead of the inverse function a(x) (it is difficult to give it in a simple form), 
we take the function x(a) as 



V(a ) - Kcl 



,l+2a 



-da — Ken 



,l+2n 



-da + Ac n 



a a 3+2n 

3(aa + 1) 



da 



2(aa+l) " "~ V J 2(aa + l) 2 
where the physical domain is V(ao) — V(a) > 0, i.e. 

-V(a ) + \Kcl f ^da + \Kc\ f ^da 



a fl 3+2n 
(aa + l) 2 



-da 







A> 



J 3(aa+l) da + 3 J (aa+l) 2 da 



where 



da 



71-1 ( -Wkrvk- 1 



m rt m— 1 



aa + 1 



a m (aa + l) ^(m-fc)!a m_fe 

da 1 

a'(cra + l) 2 a 



1 +iY-t 

a l (aa+l) [ ^ (I + 1 

\fc — 1 



In- 



On the figures|3[31 for the simplicity of presentation without loss of generality, it is assumed Cq = 1, V(a ) = 1 and 
K, a are parameters and n is chosen —1, —2, —3. 



VI. TUNNELLING IN n DECAYING COSMOLOGIES 



In the classical VSL cosmology a particle trajectory is determined through a knowledge of both a position x and 
a canonically conjugate momentum p x . In the quantum VSL cosmology the notion of trajectories loses its classical 
meaning due to the uncertainty relation [x and p x are replaced by non-commuting operators). The Wheeler-De Witt 
equation 



d 2 
ft? 



ip(a) = 



(66) 



is identical to the one-dimensional time-independent Schrddinger equation for a one half unit particle of energy E 
subject to the potential V{x) {i^(a) is known as the wave function of the universe). 

The 'particle-universe' can quantum-mcchanically tunnel through the potential barrier. Let us consider the case 
when the region beneath the barrier < a < do, is classically forbidden whereas the region a > a is classically 
allowed. 
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FIG. 3: The relationship between A and x for K = 1. The qualitative evolution of the models is represented by levels of 
A = const. The domain under the characteristic curves A(a) is classically forbidden. For n = — 1 and A > all models start 
from finite size of scale factor and expand. For n = — 1 and A > all models start from the singularity and expand. For n = —2 
and A > all models start from finite scale factor and expand. For n = — 2 and A there are oscillating universes and with a 
singularity and other universes expanding from finite scale factor. For n = —3 all models are oscillating and with a singularity. 



We can adopt a simple method of calculating an amplitude for the quantum creation of the VSL FRW universe 



from a zero size to 



_ / 3(rt+2) 



a — a y 2(„+i)A 

|VSLFRW(a )|nothing)| 2 



P = exp 



2 f x ° 



(67) 



The above formula (called Gamov's formula) gives us tunnelling probability because the quantized VSL FRW universe 
is mathematically equivalent to a one-dimensional particle of unit mass. 

As an example, let us get back to the previously considered case of the compact vacuum VSL model. The Hamil- 
tonian for this case can be obtained if we put E = into 7i.(p, x). 

The region of barrier < x < xq is classically forbidden for the zero energy particle. Therefore one can find the 
probability that particle at x — can tunnel to x — xq '■ V(xq) — 0; xq — 3 ^" +2 ^ 

After rescaling the variable x i— > x/xq = k we obtain for l|67|) 



2(n+i)A ' 



P = exp 



4 {3(n + 2)K 



h V2(n + 1)A 



(n+3)/2' 



K n/2 yjK(l - K)dK 



(68) 



where we assume K = 1, A > and — 1 < n < 0; the potential V(x) = x n+1 (l + x) has two extrema and two zeroes. 
Relation l|68|l can be rewritten to the form 



where 



F{n) 



P = exp 



3(n + 2) 
2(n+l)A 



-F(n) 



2 l(3+; 



(69) 



(70) 




X 



FIG. 4: The relationship between A and x for K — 0. The qualitative evolution of the models is represented by levels of 
A = const. The domain under the characteristic curves A(a) is classically forbidden. For n — —1 and A < all models oscillate. 
For n = — 1 and A > models evolve from a singularity to infinity. For n = — 2 and A < models oscillate. For n — — 2 and 
A > models also oscillate but without a singularity. For n — — 3 and A < models expand from a singularity to infinity. For 
n = —3 and A > all models oscillate. 



It is the most probable that the closed and vacuum VSL FRW model with — 1 < n < is born when we have 
maximum permissible energy density or least size ao- It occurs that the creation of the universe with constant c 
(n = 0) is the most probable when classical spacetime emerges via the quantum tunnelling process whereas c(a) is a 
decreasing function during the evolution of the universe. 



VII. CONCLUSIONS 



Let us assume that one takes the idea of the varying speed of light seriously as a physical effect that might have 
happened in the very early universe and today is confined to a very narrow range admissible by inaccuracy of existing 
bounds on variability of c. One of the problems arising then is to see how this modification of physics would change the 
evolution of standard Friedmann-Robertson- Walker cosmological models. So far only specific qualitative results are 
known concerning the solution of flatness and horizon problems in VSL models. In the present work we attempted to 
extend this qualitative discussion in the sense that by constructing phase space portraits of VSL cosmological models 
we were able to obtain a global view of their dynamics. In order to achieve this we have used a power law Ansatz for 
function c(t) and investigated classical Einstein equations with c allowed to be a function of time. 

The two procedure of reduction of dynamics arc proposed. In the first case we reduced the dynamics of VSL models 
to a two-dimensional Hamiltonian dynamical system with a quadratic kinetic energy form and a potential function 
depending on a generalized scale factor. In the second one we reparametrized the time variable but the scale factor 
remains a state variable. In both cases the shape of the potential and the existence of the energy integral was used 
to classify possible evolutions of VSL models. These possibilities comprise models evolving from the singularity to 
infinity, oscillatory behaviour between initial and final singularity, Einstein-de Sitter type models evolving from the 
singularity to the static world, Lemaitre-Eddington type models evolving from the static Einstein solution to infinity, 
models expanding to infinity from the finite size and finally models starting and ending with finite scale factors. 

We have dealed with the full global dynamics of VSL models. From the theoretical point of view it is important 
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FIG. 5: The relationship between A and x for K — — 1. The qualitative evolution of the models is represented by levels of 
A = const. The domain under the characteristic curves A(a) is classically forbidden. For n = —1 and A < there is no solution. 
For n = — 1 and A > models start from finite size and expand. For n = —2 there is discontinuity at x = 1 and for A < 
some models oscillate from a singularity and other models oscillate without a singularity. For n = — 2 and A > models start 
from finite size of scale factor to x = oo. For n = — 3 all models oscillate with a singularity. 

how large the class of models without horizon or with solved cosmological puzzles is. We call this class of models 
generic if its inset in the open phase is open or non-zero measure. Such a point of view is justified by the fact that if 
the solution of a cosmological puzzle is an attribute of a trajectory with a given initial conditions, it should also be 
an attribute of another trajectory which starts with neighbouring initial conditions. 

We have shown that the assumed time dependence of the speed of light leads to a uniform evolution pattern of 
VSL models on the phase space. The criteria for solving the flatness and horizon problems were formulated in terms 
of the phase space. It is an advantage of the phase space approach that one can trace the patterns of evolution for all 
possible initial conditions. We have depicted, on respective phase portraits, the regions where the flatness problem is 
solved. The models where the region of initial conditions leading to flatness and horizon problem avoidance is large 
play a distinguished role. From this perspective open (K = — 1) models with positive cosmological constant A > 
are preferred in the class of VSL FRW models filled with radiation. 

The formalism presented in this paper can be easily extended to the case where the matter content of the model is 
a mixture of different types of matter and to the case of models with shear (e.g. Bianchi I or V). 

This formalism can be also treated as a star ting point of an application of quantum cosmology to the description 
of early stages of evolution of the Universe |29ll3C| . The tunnelling rate, with an exact prefactor can be calculated to 
the first order on h for the closed VSL FRW model with a decaying variable velocity of light term c(a) The tunnelling 
probability P can be calculated in the WKB approximation given in the V ^> E limit by (|(j7fl . We consider the closed 
vacuum VSL FRW models for which potential is qualitatively classical. This implies that — 1 < n. In the interval 
— 1 < n < the probability of tunnelling increases as F(n) monotonically decreases with increasing n. It is shown 
that the highest tunnelling rate occurs for n = 0; it corresponds to the standard FRW model. 

In our work we showed the effectiveness of dynamical system methods in the investigation of VSL FRW models, 
namely 

In the distinguished class of open models with the cosmological constant the acceleration has 'transitional' character, 
i.e., there is a finite time when trajectories are in the acceleration region, and the measure of this region normalized 
to the area of phase plane is finite, even in the case of A = 0. 
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We can argue that the considered VSL models are structurally unstable (Fig.|2J because of the presence of degenerate 
critical points at infinity for n < 0. From the theoretical point of view such a situation seems to be unsatisfactory 
because in the space of all dynamics systems on the plane they form set of zero measure (the Peixoto theorem). 

The advantage of representing dynamics in terms of Hamiltonian is to discuss how trajectories with interesting 
properties are distributed on phase plane. 
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